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INTRODUCTION 



Tensor network states are powerful variational ansatze that can be used to character- 
ize the low-energy properties of quantum many-body systems on a lattice. The premise 
of tensor network approaches is to parameterize a many-body wave-function by using 
a collection of tensors connected into a network. The number of parameters required 
to specify these tensors is much smaller than the exponentially large dimension of the 
system's Hilbert space, in such a way that very large (and even infinite) lattices can 
be considered. 

Tensor network states can be broadly classified into two sets according to the geome- 
try of the underlying networks pQ. In the first set, the network reproduces the physical 
geometry of the system, as specified by the pattern of interactions in the Hamiltonian. 
For instance, the matrix product state (MPS) [2H3], an ansatz for D = 1 dimensional 
systems, consists of a collection of tensors connected into a chain; similarly, its gen- 
eralization for lattices in D > 1 dimensions, known as projected entangled pair states 
(PEPS) [MI], consists of a collection of tensors connected according to a D dimen- 
sional lattice. In contrast, a second class of tensor network states aim at reproducing 
the holographic geometry of the system. The latter spans an additional dimension used 
to parameterize the different length scales (or, equivalently, energy scales) relevant to 
the description of the many-body wave-function. Thus the multi-scale entanglement 
renormalization ansatz (MERA) [8H2"3"] for a lattice system in D dimensions consists 
of a network in D + 1 dimensions. 

The simplest and most widely employed tensor network state is the MPS. The MPS 
underlies the remarkable success of the density matrix renormalization group (DMRG) 
algorithm [2~lll27| . which for almost two decades has dominated the numerical study 
of D = 1 dimensional quantum systems, providing very accurate results for low energy 
properties. DMRG has not only become ubiquitous in condensed matter physics but 
has also found application in other fields involving quantum many-body systems, such 
as quantum chemistry [25]. Further algorithms based upon the MPS have also been 
developed, such as the time evolving block decimation (TEBD) [29, 30J algorithm, later 
reformulated as time-dependent DMRG [3T] [35] , which allows the simulation of certain 
low-energy dynamics for relatively long times. 

In this manuscript we discuss the application of the MERA to study critical systems 
in D — 1 dimensions (although most of the present formalism is directly applicable 
to D > 1 dimensions) . Given the success of MPS-based methods such as DMRG and 
TEBD for quantum systems in D = 1 dimensions, it is natural to ask whether an 
alternative approach is actually needed. A clear answer to this question is obtained by 
discussing the short-comings of the MPS representation for critical systems as well as 
by exploring the benefits of including scale invariance directly into a tensor network 
state, something that is possible with the MERA but not the MPS. 

Critical systems typically lack a characteristic length scale and are thus invariant 
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under changes of scale. One manifestation of scale invariance in critical systems is 
that correlators decay polynomially, in sharp contrast with gapped systems, where 
they decay exponentially with a characteristic correlation length. It is well-known, 
however, that an MPS with a finite bond dimension \ (where \ indicates the size of 
the tensors) can never properly capture the scale invariance of a critical state. Indeed, 
a duly optimized MPS possesses an intrinsic finite correlation length C ~ X K [33 } 134 ) . 
where k is a constant that depends on the universality class of the phase transition 
under consideration, such that correlators decay exponentially at length scales larger 
than Thus, while the MPS can accurately approximate short-range properties of 
a critical ground state, it necessarily fails to capture its properties at asymptotically 
large distances. [In practice, however, the cost of MPS-based approaches scales only 
as 0(x 3 ) with the bond dimension x- This means that one can use a very large value 
of Xi which often allows the critical behavior of a system to be accurately captured up 
to some very large length scale £.] 

On the other hand, the MERA can explicitly capture the scale invariance of criti- 
cal systems [8} fTOl \12\ [l"5Hl~7t fT9| . [22 ] . a feature that has significant advantages, both 
conceptual and practical. Tensors in the MERA are organized in layers, where each 
layer corresponds to a different length (or energy) scale. In an infinite system, scale 
invariance is then easily enforced by choosing all layers of tensors to be identical. The 
resulting ansatz is referred to as the scale-invariant MERA. Certain structural proper- 
ties of the scale-invariant MERA, such as the polynomial decay of correlators and the 
logarithmic growth of block entanglement entropy [TJ [U] , already hint at its suitability 
to represent critical systems. 

In addition, this ansatz offers direct access to the scaling operators of a critical theory, 
namely those operators that transform into themselves under scale transformations. As 
a scale- invariant/covariant object, a scaling operator must act non-trivially on a region 
of the system that has no characteristic length scale. In a (l + l)-dimensional conformal 
field theory (CFT) [55H57] (corresponding to the continuum limit of a critical quantum 
system in D = 1 spatial dimensions), the support of a scaling operator can therefore 
only be one of three possibilities: (i) an infinite line, (ii) a single point, or (iii) a 
semi-infinite line. The first type of support corresponds to a global internal symmetry 
of the CFT's Hamiltonian. The second type of support is seen to correspond to local 
scaling operators, associated to local excitations. Finally, the third type corresponds to 
non-local (or semi-local) scaling operators, associated e.g. to domain wall excitations. 
Going back to the lattice, scaling operators are distorted by the presence of a finite 
lattice spacing, but they can still be directly extracted from the scale-invariant MERA. 
Thus, on the lattice, (i) a global internal symmetry is implemented by an infinite 
string of identical single-site operators; (ii) local scaling operators are now supported 
on a small number of sites (the specific number depends on the MERA scheme); and 
(iii) non-local operators mix elements of the two previous objects: they consist of a 
semi-infinite string of identical single-site operators (the same ones that implement 
an internal symmetry) completed with some local operator at the end of the string. 
Importantly, the scaling dimensions and fusion rules of the scaling operators on the 
lattice are the same as in the continuum. As a result, a relatively simple and inexpensive 
calculation with the scale-invariant MERA can be used to obtain remarkably accurate 
estimates of the conformal data characterizing the underlying CFT. 

The rest of the manuscript is organized in sections as follows. Sect. |H| introduces 
the key aspects of entanglement rcnormalization, which is the renormalization group 
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transformation for lattice models on which the MERA is based. We describe how local 
operators transform under the coarse-graining transformation and discuss basic aspects 
of the MERA, including the evaluation of the expectation value of local observables, 
and briefly compare different MERA implementations in D — 1 dimensions. Sect. |III| 
addresses the role of spatial and internal symmetries in tensor network states, and 
compares how different symmetries can be enforced in MPS/PEPS and MERA. While 
translation invariance is naturally implemented in an MPS and PEPS, it can only be 
approximately enforced on the MERA. Sect. |IV| specializes on the implementation 
of scale invariance in the MERA and discusses how the scaling operators of a criti- 
cal theory can be extracted from it. In Sect. [V]wc demonstrate the performance of 
the scale-invariant MERA for a number of critical quantum spin chains, namely Ising, 
Potts, quantum XX and Heisenberg models. Specifically, Sect. |V A| compares ground 
state energy and two-point correlators obtained with MERA and MPS. Interestingly, 
MPS and MERA approaches seem to complement each other. For a given computa- 
tional cost, the MPS is seen to provide more accurate estimates of the expectation 
value of a local observable, such as the ground state energy. However, the MERA is 
seen to provide a better characterization of long-range properties, such as two-point 
correlators at long distances. The advantages of the scale-invariant MERA are then 
further illustrated in Sect. |VB| by extracting, in the concrete context of the quantum 
Ising model, the conformal data (including scaling dimensions, fusion rules for scaling 
operators and central charge) of the underlying CFT. 




FIG. 1. (i) The coarse-graining transformation U, a specific implementation of entanglement 
renormalization, is comprised of isometries w and disentanglers u and maps blocks of three 
sites from the initial lattice C into a site of the coarser lattice £.' . (ii) The tensors w and u 
are constrained to be isometric, see also Eq. [3] 
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II. ENTANGLEMENT RENORMALIZATION AND THE MERA 

In this section we first recall a few basic aspects of entanglement renormalization 
(ER), the coarse-graining transformation upon which the MERA is based. We then 
introduce the MERA and review a few of its features. 



A. Foundations of Entanglement Renormalization 

For concreteness, we mostly consider a specific implementation of ER that produces 
the so-called ternary MERA (where three sites are coarse-grained into one effective 



site). In Sect. II C we also discuss other MERA schemes. 

Let C denote a D = 1 dimensional lattice made of N sites, where each site is described 
by a Hilbert space V of finite dimension d, so that the vector space of the lattice is 
Yc — V® N . We consider a coarse-graining transformation U that maps lattice £ to a 
coarser lattice £ 

U^-.Vc^Yo. (1) 

where £ is made of N/3 sites, each with a vector space V of dimension Xi so that 
V.C' = Y'® N , and where transformation U decomposes into local transformations, 
known as disentanglers u and isometries w, 



u 



according to Fig. [iji). More specifically, if we partition the initial lattice C into 
blocks of three sites, then the disentanglers u are first applied across the boundaries of 
neighboring blocks, followed by the isometries w, which map each block of three sites 
into a single effective site of the coarser lattice £ . Disentanglers and isometries are 
required to satisfy isometric constraints, namely 



,t„, - ir'®2 „„t„„ - it' (3) 



U ' U = 1 , w'w — 



where I and F are the identity operator in V and V', respectively, see Fig. |T|(ii) . 
Note that, by construction, the disentangler u is also unitary, that is uu^ = I® 2 . The 
dimension \ of the Hilbert space V can be chosen to be different than d, provided that 
X < d 3 (as demanded by the above isometric constraint on w). In general, choosing 
a larger dimension i.e. retaining a larger effective Hilbert space V for each coarse- 
grained site, yields a more accurate RG transformation, one that better preserves the 
low energy properties of the system. 

An important property of the coarse-graining transformation U is that it preserves 
locality. Let 0(7*1, be a local operator defined on two contiguous sites {r\,r2) of 
lattice £ Under coarse-graining, the operator o(ri,r 2 ) becomes 

U^o{r 1} r 2 )U = ■ ■ ■ ® I' ® o'Cr'n^a) ® I' ® ■ ■ ■ , (4) 

where the only non-trivial part of lPo(ri,r2)U is a new operator o'ir^r^) supported 
on two contiguous sites (r^r^) of lattice £ . Notice that operator o'{r'i,r'<^) remains 
local (i.e., it is supported on two sites) thanks to both the specific decomposition of U 
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• • • • • • • • • -ca- • • L 

o(r„r 2 )"" r > r 2 ; 3 r * 




FIG. 2. (i) Under coarse-graining with entanglement renormalization, an operator o(r\,r2) 
supported on two sites of lattice C is transformed into a new operator o'(ri,r' 2 ) supported 
on two sites of the coarser lattice see Eq. [4] The coarse-graining of local operators can 
be implemented directly via the (iia) left, (iib) center and (iic) right ascending superopera- 
tors, denoted Al,Ac and Ar respectively. Notice that the coarse-graining of o(ri,r2) in (i) 
corresponds to application of the left ascending superoperator Al- 



into disentanglers and isometries, and the isometric constraints of these tensors, Eq. |3j 
which ensure that most of the tensors of U in U'o(ri, T2)U annihilate to identity with 
their conjugates in U\ as shown Fig. [2|i). In view of this fact, it is most convenient to 
introduce left, center and right ascending superoperators {„4i, .4c", „4#}, as shown Fig. 
[2^ii), which directly produce the two-site coarse-grained operator d from the two-site 
operator o, as given by one of the following, 

d (r[,r' 2 ) = A L (ofa.ra)), 
d (r£,r£) = Ac (o(r 2 ,r 3 )), 

d(r' 1 ,r' 2 )=A R (o(r 3 ,r 4 )), (5) 



where the specific choice of ascending superoperator to be used depends on the location 
of the operator o on the lattice C. 

We may now concatenate the coarse-graining transformation a number T of times 
to obtain a sequence of coarser lattices, 

£ [o]^ £ [i]^...^ 11 / ;[T] ) (6) 

where we use superscripts in square brackets to denote the level of coarse-graining, 
with the initial lattice = C. Then, for any local (i.e., two-site) operator = o 
defined on u-°\ the transformations {J/' r '} generate a sequence of local coarse-grained 
operators {oM}, defined on lattices {£' r '}, 




FIG. 3. The ternary MERA for a lattice C [0] of N = 54 sites. Each layer U [t] of the MERA 
can be interpreted as a coarse-graining transformation between an initial lattice £' r ' and a 
coarser lattice £' T+1 '. The past causal cone of two sites (r\,r2) in lattice £' ' is shaded. 



B. Foundations of the MERA 



We have just seen that the ER transformation U can be used to coarse-grain local 
operators, producing a renormalization group (RG) flow for local operators, Eq. [7] 
As a linear (isometric) map from Yc to Wc, U can of course also be used to coarse- 
grain quantum states. More important for us, however, is to consider an inverse 
RG flow of states. Let us assume that we have the sequence of ER transformations 
{U [0] ,U [1 \. . . ,U [T - 1] } which act on an initial lattice C [0] of N sites to eventually 
produce coarse-grained lattice £l T h Then for a quantum state \ip^) defined on lattice 
C\ T \ the transformation C/P" -1 ! can be used to obtain a new state |-0' T ~ 1 '), 



^[T-l]\ = u[T -X] 



1> 



[T] 



(8) 
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defined on the finer lattice C^ T ^. Through iteration of Eq. |8l one can further obtain 
increasingly fine-grained states, eventually reaching a state defined on £)°\ 



^ [0] ) -c/ [01 c/ [11 ---t/ [T - 11 i> 



(9) 



Let us assume that the number of levels T is chosen T « log 3 (A), such that the 
maximally coarse-grained lattice contains a small number of sites and hence the 
state can also be described with a small number of parameters. Then the multi- 

scale entanglement renormalization ansatz (MERA) is the class of states that 
can be represented as Eq. |9j for some choice of {t/' ', C/W, . . . , U^ T ^^} and |"0' T ')- For 
instance, Fig. [3] depicts the MERA, organized into T = 3 layers, for a state on 
a lattice of N = 54 sites, 



^[o]\ = u [o] u m u m ^[3] 



(10) 



Let us now count variational parameters. We will assume for simplicity that for any 
value of r = 0, 1, • ■ • , T— 1, the dimension of the vector space is \. Recall that the 
transformations {C/M} themselves are comprised of local tensors, the disentanglers u 
and isometries w, each specified by \ A parameters. Since in an A^-site lattice there are 
O(N) disentanglers u and isometrics w (d istri buted in (9(log(A)) layers), the MERA 
depends on 0(Nx 4 ) parameters. In Sect. 
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we will see that, through incorporation 
of spatial symmetries into the ansatz, this number can be reduced to 0(x 4 ), which is 
independent of N and allows for the study of infinite systems. 

We have therefore established that the MERA can be specified with a number of 
parameters that is much smaller than the dimension of the Hilbert space Yc, which 
grows exponentially in the number N of sites. However, for this ansatz to be useful, 
we also need to be able to efficiently extract information about in Eq. [9] For any 

local operator (?"!, r 2 ), and due to the very peculiar causal structure of its tensor 
network, it is actually possible to efficiently compute the expectation value 



-,[°] 



(n,ri) 



[0] 



(11) 



from the MERA. Let us define the past causal cone of a site in lattice C 1 ' 7 ' as the set 
of tensors and indices that can affect the state on that site. By construction, in a 
MERA the causal cone of any site of is seen to involve just a constant (that is, 
independent of N) number of sites of any other lattice C} T 1 for r' > t, a property that 
we refer to by saying that the past causal cone has bounded 'width'. Fig. [3] displays 
the past causal cone of two sites (ri,r 2 ) in a ternary MERA, which only involves two 
sites of every lattice C^. This property allows for the efficient computation of local 
reduced density matrix /^(ri,^), from which the expectation value 



J ! 



(n,n))=tr (o^Crx.nJpMCn.ra)) 



(12) 



can be obtained. The computation of local reduced density matrices is simpli- 
fied through the introduction of left, center and right descending superoperators, 
{1>l, T> c, T> r}, which are the adjoints of the ascending superoperators of Eq. [5] Let us 
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(i) 



(iia) 




r i r 2 7 3 '4 



(iib) 



(iic) 



P 1 ' l >(r 2 ,r 3 )-V c (p"(r;,r;)) 



p lt - n (r„r 4 ) = V M (p^(r;,r^) 





FIG. 4. (i) The causal cone (shaded) of four sites (ri,ra,r3,r4) in lattice £' T_1 ' involves two 
sites (ri,r' 2 ) in lattice C^ T \ Starting from the reduced density matrix p' T '(»"i, 7-2) on lattice 
£' T ' , the reduced density matrix on any pair of contiguous from (n , r 2 , r$ , n) can be obtained 
using the (iia) left, (iib) center and (iic) right descending superoperators, denoted T>l,T>c,'Dr 
respectively. 



assume that we have the density matrix (r[, r 2 ) describing the state on two contigu- 
ous sites (r[, r 2 ) of lattice . Then, as shown Fig. [4J the descending superoperators 
may be used to compute the two-site reduced density matrix p^ T_1 l on certain sites of 
the lattice C^" 1 ^, 

P [T - 1] (n,r 2 )=V L (pMOVa)), 
P lT - 1] (r 2 ,r 3 )=V c (pM(riy 2 )), 

P [T - 1] (rs,r i )=V R (pW(r' i y 2 )), (13) 

Thus, through repeated use of the appropriate descending superoperator in Eq. [T3j we 
can compute the reduced density matrix (n , r 2 ) of any two contiguous sites (n , r 2 ) 
of the original lattice by 'lowering' the density matrix through the appropriate 
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causal cone. For instance, the reduced density matrix p[°' (rj., r%) on the two sites 
(ri,r 2 ) of lattice shown in Fig. [3] can be computed as 



p [0] (r 1 ,r 2 )=V c (V L {V. 



1> 



(3) 



(3) 



(14) 
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The cost of calculating a local reduced density matrix , as in the example of Eq 
is proportional to the number T w log(A) of layers in the MERA, hence scales with 
system size N as 0(\og(N)). Two point correlators can also be evaluated using similar 
manipulations, see Ref. 1151 



C. Choice of MERA Scheme 



There are many possible ways of implementing the MERA in D = 1 dimensions, 
of which the ternary MERA described in Sect. |IIB| is just one example. Fig. [5] 
displays a ternary MERA together with two other possible implementations: the binary 
MERA scheme (in terms of which the first ER proposals [U [9] were formulated) and a 
modified binary MERA scheme (with half the amount of disentanglers as the previous 
binary scheme). While all MERA schemes function similarly on a conceptual level, 
the computational efficiency may differ between schemes. For instance, the cost of 
evaluating the expectation value of local observables as described Sect. |IIB| scales, 
in terms of the bond dimension X: as 0(x 9 ) f° r the binary MERA, as 0(x 8 ) for the 
ternary MERA and as 0(x 7 ) for the modified binary MERA. On the other hand, 
the binary MERA scheme has more disentangling power than either the ternary or 
modified binary schemes and, for any given x, produces a more accurate representation 
of ground states. It is therefore not obvious which MERA implementation will give the 
best numeric results for a fixed computational budget. However, a direct comparison 
of performance, see Sect. |Bj shows that the modified binary scheme of Fig. [5|ni) is the 
most efficient scheme. Consequently, this scheme is used for the obtaining the numeric 
results presented in Sect. [V] However, for the sake of simplicity, we shall continue to 
discuss theoretical aspects of MERA in terms of the ternary MERA. 



III. SYMMETRIES IN TENSOR NETWORK STATES 



Consider a many-body state that is invariant under some symmetry transformation. 
In approximating this state with a tensor network state, we would like to preserve 
the original symmetry as much as possible. In this section we examine the types of 
symmetries that can be enforced upon the MPS and PEPS, and upon the MERA. 
We also examine whether the presence of the symmetry can be exploited for compu- 
tational gain. We begin by discussing spatial symmetries, followed by global internal 
symmetries. The results are summarized in Table [IJ 



A. Spatial Symmetries 



For simplicity we discuss only two typical spatial symmetries: invariance under trans- 
lations in homogeneous systems and invariance under changes of scale in e.g. critical 
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TABLE I. Several symmetries that can be exactly enforced and/or whose presence can be 
exploited for computational gain with MPS/PEPS and MERA algorithms. 



Symmetries with MPS/PEPS: 





Enforceable 


Exploitable 


Translation invariance 
Scale invariance 
Internal symmetries 
(e.g. Za,U(l),SU(2)) 


Yes 
Unlikely 

Yes 


Yes, cost: O(N) -> O(l) 
Unlikely 

Yes 


Symmetries with MERA: 




Enforceable 


Exploitable 


Translation invariance 
Scale invariance 
Internal symmetries 
(e.g. Z 2 ,U{1),SU{2)) 


Unknown 
Yes 

Yes 


Yes, cost: O(N) -> OQog(N)) 
Yes, cost: 0(log(AT)) -> O(l) 

Yes 



systems. Let us first consider them in the context of the MPS for D = 1 dimensions 
and PEPS for D > 1 dimensions. 

In an inhomogeneous MPS/PEPS, one associates a different tensor A? to each site r 
of the lattice C. Hence, in a lattice made of N sites, the total number of tensors in the 
tensor network is also N . In the presence of translation invariance (either in a finite 
system with periodic boundary conditions or in an infinite system), this symmetry can 
be incorporated into the MPS/PEPS by choosing all the tensors to be a copy of the 
same tensor A, i.e. A? = A, see Fig. [6|l) for an MPS. Translation invariance is in 
this way exactly preserved. It can also be exploited to reduce computational cost from 
O(N) to 0(1). 

On the other hand, it is not clear how scale invariance could be enforced in these 
tensor networks. For an MPS this is unlikely to be possible at all because, as we 
mentioned earlier, a finite bond dimension \ already implies the presence of an effective 
finite correlation length C ~ X K 031 IS] . 

Let us now consider spatial symmetries in the MERA. Recall that a generic MERA 
on an N site lattice is arranged into T s» log(iV) layers of tensors, and contains 0(N) 
different tensors, as depicted in Fig. [6]jiia). Suppose now that the state to be approx- 
imated by the MERA is translation invariant. Then we can choose all the tensors in 
each layer to be the same, so that layer is characterized by a single pair of ten- 
sors u^ T l and w^- T \ see Fig. [6|iib). In this way translation invariance can be exploited 
to reduce computational costs from 0{N) to 0(log(AT)). Notice, however, that this 
choice of tensors does not enforce translation invariance, because the structure of the 
coarse-graining is not homogeneous (different sites are positioned in incquivalent posi- 
tions with respect to the disentanglers and isometries). The final effect is examined in 
Fig. [7J A MERA characterized by a single pair of tensors u 1 ^ and for each layer, 
where these tensors are filled with random coefficients (compatible with the isometric 
constraints of Eq. [3]), is highly non-translation invariant, with e.g. oscillations in the 
expectation value of the energy of the order of 0.1 for the Hamiltonian i?i s i ng of Eq. 
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>i ' 2 

FIG. 5. Three different MERA schemes for a D = 1 dimensional lattice. An example of a 
causal cone is shaded for each scheme, (i) The binary MERA scheme, based upon a 2-to-l 
coarse-graining step, has a causal width of three sites and a cost of contraction that scales 
with the bond dimension x as 0(x 9 )- (h) The ternary MERA scheme, based upon a 3-to-l 
coarse-graining step, has a causal width of two sites and a cost of contraction that scales as 
0(x 8 )- (iii) The modified binary MERA scheme, equivalent to the binary MERA scheme with 
every second disentangler removed, has a causal width of two sites and a cost of contraction 
that scales as 0(x 7 )- 
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(ia) (ib) A 

A A A A A A A ^ ^ \. 




FIG. 6. (ia) An inhomogeneous MPS has an independent tensor Ap associated to each lattice 
site, (ib) All the tensors in a translation invariant MPS are chosen to be copies of the same 
tensor A. (iia) A generic MERA on an N site lattice contains O(N) different tensors, (iib) 
Translation invariance can be exploited by choosing the tensors in each layer C/' T ' of the MERA 
as copies of a single unique disentangler u' r ' and isometry w' r '. For an N site lattice, this 
MERA contains 0(log N) different tensors, (iic) Scale invariance can be incorporated into the 
MERA by further enforcing all layers to be identical, hence the entire MERA is described by 
a single u and w. 



[36] Still, these violations of translation invariance decrease significantly once the ten- 
sors are optimized so as to minimize the expectation value of the translation invariant 
Hamiltonian -f/ising- Indeed, they become of order 10 -5 for \ — 4 and decrease with 
increasing \- P n practice one can efficiently average the expectation value of a local 
observable over all possible lattice positions in order to further reduce the effect of 
these small violations of translation invariance]. We conclude that translation invari- 
ance can be exploited to reduce computational costs, but it can only be reproduced 
approximately. It is not known whether it can be enforced exactly. 

Instead, enforcing scale invariance in the MERA is straightforward. This is accom- 
plished by choosing all disentanglers and isometries to be copies of a single pair u and 
w, see Fig. [6]jiic), which further reduces the number of parameters and the computa- 
tional cost of MERA algorithms from 0(\og(N)) to O(l), allowing infinite systems to 
be considered. The scale-invariant MERA will be discussed in more detail Sect. HVl 

To summarize, in the MPS/PEPS we can enforce and exploit translation invariance 
but not scale invariance, whereas in the MERA we can enforce and exploit scale in- 
variance but only exploit (i.e., we cannot enforce) translation invariance. Thus both 
MPS/PEPS and MERA have potential advantages over each other, depending on 
whether exact translation invariance or exact scale invariance is more important for 
the problem under consideration. 
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FIG. 7. We investigate translation invariance in the scale-invariant MERA by comparing the 
bond energy E(r) over 30 contiguous lattice sites with the average bond energy E from all 
sites, as measured with the critical Ising Hamiltonian, -Hi s i ng , of Eq. |36| For a randomly 
initialized X = 4 scale-invariant MERA, the large fluctuations of bond energies indicate the 
state is highly non-translationally invariant. Once the MERA has been optimized for the 
ground state of //i s i ng , it more closely approximates translation invariance; bond energies now 
differ from the average by less than 0.1%. As the bond dimension \ of the MERA is increased, 
the optimized wavefunction better approximates translation invariance; for \ = 16 the bond 
energies differ from the average energy by less than 0.001%. 



B. Global Internal Symmetries 



A second important class of symmetries are those involving internal degrees of free- 
dom, such as Z2 spin flips and U{1) or SU(2) spin rotations simultaneously applied 
on all the sites of a spin model. Such symmetries can be enforced and exploited in all 
tensor networks. 

Let us assume that the Hamiltonian H of our lattice model is invariant under a 
symmetry group Q, 

T g HT\ = H, Vg e G, (15) 

where T g = ■ ■ ■ V g ® V g ® V g ■ ■ ■ is an infinite string of copies of a matrix V g , with V g a 
unitary representation of Q, and let \ip) be the ground state of H, which we will assume 
to have the same symmetry, i.e. T g \tp) = \ip) (or, more generally, T g \tp) = e l ^\4>))- 
We can then ensure that the symmetry is also exactly preserved in a tensor network 
approximation to \ip) by using symmetry preserving tensors 38, 39]. For instance, for 



1G 




FIG. 8. In order to preserve a (global) symmetry specified by symmetry group Q, the tensors 
u and w comprising the MERA are chosen to be invariant under the action of a unitary 
representation V g of symmetry group Q, see also Eq. |16| 

the MERA, we choose the disentanglers u and isometries w such that, 

{V g ® V g ) u {V g ® V g y = u, 

{Vg ® Vg ® Vj) Itf (V^ = (16) 

where acting on different indices may actually denote different (in general, reducible) 
representations of G, see also Fig. [8] The use of symmetry preserving tensors implies 
that the tensors are block diagonal when viewed in a certain basis and thus contain 
less free parameters than generic tensors. This reduction in the number of parameters 
can be exploited to significantly decrease computational costs. Symmetries, and in 
particular a truncated version of the operator r g , also play an important role in the 
description of non-local scaling operators, as discussed in Sect. |IVD| 

IV. SCALE-INVARIANT MERA 

We have already introduced the scale-invariant MERA: in a lattice C with an infinite 
number of sites, N — > oo, it consists of infinitely many layers of tensors, where all the 
disentanglers and isometries are copies of a unique pair u and w. In this section we 
enumerate two significant structural properties of the scale-invariant MERA and review 
how one can compute a local reduced density matrix, from which the expectation value 
of a local operator can be evaluated. Then we discuss the three types of scale-invariant 
(or covariant) objects one can extract from it. 

A. Basic Properties 

Two basic features of the scale-invariant MERA in D = 1 dimensions match well- 
known properties of the ground state of a critical system. Firstly, the entanglement 
entropy Sl of a block of L contiguous sites can be seen to scale as the logarithm of L 
[9] , which is compatible with the critical scaling [40l E] , 



S L w|log(£), 



(17) 
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FIG. 9. A scale-invariant MERA consists of a number M of transitional layers, here M — 2 
transitional layers f/' ' and J/' 1 ', followed by an infinite number identical of scale-invariant 
layers U. Recall that each layer C/' T ' is comprised of local isometric tensors, the disentanglers 
and isometries ui' T ', as depicted Fig. |l[i). 



where c is the central charge of the CFT. Secondly, correlation functions can be seen 
to decay polynomially [S], 

(o(n)o(r 2 )) « 1 (18) 
|n - r 2 \ q 

as it is also expected of critical correlators. Interestingly, these two properties of the 
scale-invariant MERA follow from simple geometric considerations, namely by studying 
minimally connected regions and geodesic paths in the (discrete) holographic geometry 
generated by the tensor network pQ. 



B. Transitional Layers 

In a practical computation (see Sect. |A| it is customary to consider a scale- invariant 
MERA with some small number M of translational layers {U^°\ • • • , [/[ M_1 l}, which 
are characterized by M pairs of tensors {(u^\w^), • • • , (ut M-1 J , u>[ M_1 J)} that are 
chosen independently of the single pair (it, w) characterizing the rest of layers (see Fig. 
Mi) for an example with M = 2). These transitional layers serve two main purposes. 
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Firstly, they allow one to choose the bond dimension \ °f the scale-invariant layers 
independent of the local dimension d of the sites in the orig inal lattice £[°1 . 

Secondly they also allow to diminish the effect of RG irrelevant terms in the critical 
Hamiltonian H of the system. Such terms violate scale invariance but become less and 
less important at larger length scales. The number M of transitional layers required 
depends on the amplitude and scaling dimensions of the irrelevant operators present 
in H, and is often determined by trial and error. For the sake of simplicity, in the rest 
of this section we shall focus on the case of a purely scale-invariant MERA with no 
transitional layers. 



C. Local Density Matrix 



The computation of the (average) local density matrix 



0x> 

\ r=l 



p^lirnj , \ ) (IT) 

for two contiguous sites of lattice £ is of central importance in the present formalism. 
The density matrix p is required both in order to extract the expectation value of a 
local operator o from the scale-invariant MERA and to optimize its tensors so as to 
approximate the ground state of a critical Hamiltonian H . 

Here we consider the evaluation of the expectation value (o(r,r + 1)) of a local 



observable o(r, r + 1). As discussed in Sect. Ill A the MERA in not manifestly trans- 
lation invariant. Thus the expectation value (o(r, r +1)) can artificially vary with the 
position of site r in the lattice. To mitigate this effect, rather than evaluating the 
expectation value at a particular lattice position r we will instead evaluate the average 
expectation value over all lattice sites, 



(o)= lim _y>( r , r + l)) . (20) 




Notice that this average expectation value can be expressed in terms of the average 
two-site reduced density matrix p introduced in Eq. |19[ 

(o) =tr( o(r,r + l) p). (21) 

In Sect. |II B| we described the use of the left, center and right descending superoper- 
ators, {T>l,T>c, ~Dr}, to compute the reduced density matrix pl°l (r, r + 1) from a finite 
MERA. In particular, it was argued that obtaining the density matrix p>- '(r,r + 1) 
required application of a specific sequence of left, center and right descending super- 



operators that depended on the causal cone associated to sites (r,r + 1), see Eq. 14 
for an example. The average density matrix p can be seen to follow from using the 
average descending superoperator V, defined as 



V = i (V L + V c + V R ) , 



(22) 
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in order to descend the density matrix through the 'average' causal cone. That is, 
given the average density matrix pM at level r, the average density matrix pt r_1 l at 
lower level t — 1 is obtained as 

p^- 1 ! = V (p M ) . (23) 

In an infinite system, n — ¥ oo, the MERA has T —¥ oo layers, and the average density 
matrix p is obtained from 

p=hjpo...oP)(^), (24) 

T times 

where p is simply the dominant eigenoperator of the descending superoperator T>, which 
is independent of pI T l . As a manifestation of scale invariance, this is also the two-site 
density matrix of any coarse-grained lattice C^ T \ that is pM = p for any r > 0. More 
details on the computation of p can be found in Sect. [X] 



D. Scale-invariant Objects 



The scale-invariant MERA offers direct access to objects of a critical theory that are 
invariant (more generally, covariant) under a change of scale. In D = 1 dimensions we 
can identify three classes of such objects, as depicted Fig. 
them in some detail. 



10 'i-iii). Next we discuss 



1. Symmetry Transformations 

Let us assume that the critical ground state is invariant under some symmetry group 
Q, as implemented by the symmetry transformations T g = • • • V g <S)V g (S)V g ■ ■ ■ introduced 
in Eq. |15[ where g € Q . The first type of objects that transform into themselves under 
changes of scale correspond precisely to the infinite strings T g , for which we have 

r 9 ^> r g (25) 

Indeed, if we use symmetry preserving tensors in the MERA, as per Eq. [THJ then 
the string T g commutes with each layer U of the MERA, T g U — UT g , or equivalently 
the string T g remains invariant under coarse-graining, U^T g U = T g , as shown in Fig. 
Mi). 



2. Local Scaling Operators 

The second class of objects with a simple transformation rule under coarse-graining, 
which can be easily extracted from the scale-invariant MERA, are local scaling opera- 
tors <p a , fulfilling 

4>a > K 4>a, (26) 
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(iv) 



o' = S(o) 



\ 



5 



-o 



(v) 



°' = s g {°) 




FIG. 10. The three classes of scale-invariant/covariant objects for a D = 1 dimensional lattice 
include (i) an infinite string of V a , with V B a unitary representation of symmetry group Q of 
the system, (ii) local operators and (iii) non-local operators, which consist of a local operator 
with a semi-infinite 'tail' of V g . (iv) Scaling superoperator S for local operators and (v) scaling 
superoperator S g for non-local operators. 



where X a is some constant. 

For simplicity, here we focus on one-site scaling operators [One can also compute 
two-site scaling operators, but they lead to the same scaling dimensions and fusion 
rules]. As depicted in Fig. 10 n), a one-site operator o located at certain points on 



the lattice is coarse-grained into another one-site operator d . This coarse-graining can 
directly be implemented with the one-site ascending superoperator, which we call the 
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(one-site) scaling superoperator S in the scale-invariant setting, 

o'=S(o), (27) 
see also Fig. [lOjiv). Iteration produces an RG flow for one-site operators: 

o A d A o" ■■■ . (28) 
The scaling operators 4> a and their corresponding scaling dimensions A„, 

S(4> a ) = X a ^ a , A Q = -log 3 A Q , (29) 
can then be obtained by simply diagonalizing the scaling superoperator S |141 116j . 

3. Non-Local Scaling Operators 

Let us assume again that the critical ground state represented by the scale-invariant 
MERA is invariant under a symmetry group G, as implemented by the symmetry 



transformations r g = • • • V g <S> V g <£> V g ■ ■ ■ introduced in Eq. 15 where g e Q, and that 
the tensors of the MERA have been chosen to preserve this symmetry, as per Eq. [16| 
We can then identify a third class of objects with a simple transformation rule under 
changes of scale, namely non-local scaling operators <\f g a , to be defined below, which 
fulfill 



' 1 A,.... o;;.„. (30) 



where A SiQ is some constant. 

To see how these scaling operators come about [32], let us first introduce non-local 
operators o g of the form, 

^r> , rt = — v B ®v g 2>v g (3i) 



where T g is a semi-infinite string made of copies of V g and o is a one-site operator 
attached to the open end of 1^. Notice that, under coarse-graining, can be mapped 
into another non-local operator o g ' of the same type, 

o2 = rs® ^ < = r<® ', (32) 

since the semi-infinite string T g commutes with the coarse-graining everywhere except 
at its open end, as illustrated in Fig. [ToViii) . Thus we can study the sequence of 



coarse-grained non-local operators o g — > o g ' — > o g " ■ • • by just coarse-graining the 
local operator o with the modified one-site scaling superoperator S g of Fig. 10 V), 



o — > o — > o ■ ■ ■ . (33) 
In particular we can diagonalize the modified scaling superoperator S g , 

Sg(4>g,a) = ^g, a 4>g,a , &g,a = - fog 3 \,a , (34) 
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to obtain non-local scaling operators tf>g a of the form 

C = r«®<^ )Q . (35) 

Notice that for g — I we recover the local scaling operators (j> a of Eq. [29] 

Importantly, the scaling dimensions Aq, and A g Q , of the both local and non-local 
scaling operators </> a and <f>g a (as well as their operator product expansion coefficients, 
see Refll6j) are the same in the lattice than in the continuum. Therefore by extracting 
properties of the scaling operators on the lattice, we can characterize the CFT that 
describes the critical theory in the continuum. As demonstrated by the benchmark 



results of Sect. VB[ a relatively simple and inexpensive MERA simulation can actually 



be used to obtain remarkably accurate conformal data of the underlying CFT. 



V. BENCHMARK RESULTS 

In this section we benchmark the performance of the scale-invariant MERA by ap- 
plying it to study of the ground state of several well-known quantum spin chains at 
criticality. The models we analyze are the critical Ising model [42] [43] , the critical 
three-state Potts models [H], the XX model [J5] and a Heisenberg zig-zag chain [JS] 
(the Heisenberg model with a next-nearest neighbor coupling), corresponding to the 
following Hamiltonians: 



-^Ising 


=£ 

r 


(z{r)-X(r)X{r + l)) 




(36) 


-ffpotts 


= £ 

r 


' Z{r) - X(r)X\r + 1) - A f (r 


)X(r + l)) 


(37) 


Hxx 


=£ 
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r X(r)X(r + l) + Y(r)Y(r + l) 
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(38) 
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1-2)) 


(39) 



/-10 






fo 


1 




2 


•) 


, x = 








i) 


\ 






K 1 








where X, Y, Z are Pauli matrices, S = [X, Y, Z], and where Z, X are three-state Potts 
spin matrices given by 



(40) 



The next-nearest neighbor coupling in the Heisenberg zig-zag chain is set at the critical 
value Ji = 0.24116; at this value the model is scale-invariant 0B]. Note that, although 
the standard Heisenberg model (with J2 = 0) is quantum critical, the Hamiltonian 
contains a marginally irrelevant contribution that breaks scale invariance, which here 
we remove by adding the next-nearest neighbor coupling. 

In the present calculation, we have used the modified binary MERA scheme, de- 
picted in Fig. , with either M — 2 or 3 transition layers. This ansatz is optimized 
by minimizing the expectation value of the energy density for each of the above Hamil- 
tonians, by using the optimization algorithm described in Sect. [A] In terms of the 
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bond dimension Xmera > the cost of optimizing the modified binary MERA scales as 
0(Xm ERA ), but can be reduced to O(Xmera) through use of an approximation in the 
tensor network contractions as described Sect. [C] We have computed the ground states 
of the four models over a range of values of Xmera up to Xmera — 150. Each simulation 
took under a week on a 3GHz dual-core workstation with 32Gb of RAM. 

For purposes of comparison, we have also computed the ground state of the four 
critical spin chains using an infinite, translation invariant (with a 4-site unit cell) MPS. 
The MPS tensors are optimized with a variational approach similar to the iDMRG 
algorithm 47J. The computational cost scales with the bond dimension x M ps of the 
MPS as O(Xmps)- We have computed the ground states of the four models over a range 
of values of Xmps up to Xmps = 1536. 

In both the MERA and MPS ca lculations we have employed symmetry preserving 
tensors, as described in Sect. IIIB to enforce (some of) the global internal symmetries 
of these models. Specifically Z2 symmetric tensors have been used for the Ising model; 
Z 3 symmetric tensors have been used for the Potts model (Z3 is a subgroup of the 
full S3 symmetry of this model); and £7(1) symmetric tensors have been used for both 
the quantum XX and Heisenberg zig-zag chains (again, U(l) is a subgroup of the full 
SU(2) symmetry of the Heisenberg zig-zag chain). 

In the first part of the benchmark, Sect. |V A| we compare ground energy and two- 
point correlators obtained from MERA and MPS, and discuss the relative merits of 
each approach. Then in Sect. |VB| we demonstrate the extraction of conformal data 
from the scale-invariant MERA for the critical Ising model, following the approaches 
of Refs.QHHl]. 



A. Comparison with MPS 

Here we compare the performances of MPS and scale-invariant MERA for the com- 
putation of ground state energy and two-point correlators. 



1. Ground Energy 

For both Ising and quantum XX models, the exact ground energy per site is E — 
— 4/7r, while for the three-state Potts and Heisenberg zig-zag chains we use an MPS 
with Xmps = 1536 to estimate the ground energy per site at -Ep tts = —2.4359911239(1) 
and ^Hcis.z.z. ps -1.607784749(1). 

Fig. ITT] displays the relative error in the ground state energy per site, AE = (E ex!ict — 
£? numer ic)/-Eexact for the models under consideration over a range of bond dimensions 
X, for both MERA and MPS. This figure reveals a number of interesting similarities 
between results obtained with MERA and MPS. Recall that the central charge c for 
these models is 

1 4 

Clsing = ^' c Potts = C XX = 1; c Hcis.Z.Z. = 1, (41) 

Then a first observation is that for both MERA and MPS, for a given bond dimension \ 
the larger the central charge c the larger the error in the energy is. A second similarity 



24 



1 o 



4 . 



< 
o 



10" 



1 0" 



10 



lu 10 
>■. 

01 

,n ID" 2 



01 



1 0" 



1 0" 



1 0" 



- 1 



10 



(i)MERA 



fc- ^: 



• 


Is ing model 


+ 


3-state Potts model 


X 


XX model 





Heisenberg Zig-Zag 



.....^ 



10 20 40 

Bond Dimension, Xmera 



200 



(ii)MPS 



... 



1 



1 



Bond Dimension, Xmps 



FIG. 11. Relative energy error AE in the ground state of critical spin models as a function of 
the tensor network bond dimension \y comparing (i) the scale-invariant MERA with (ii) an 
infinite MPS. In all cases the energy error AE appears to scale polynomially in \ in accordance 
with Eq. 1421 



is that for both MERA and MPS the error AE in the energy scales polynomially in x, 
i.e. to a good approximation, 

AE = ax~ b . (42) 

A linear fit in Fig. [IT] produced the estimates for the coefficients a and b displayed 
in Table [TTJ In the large x regime, the error AE is dominated by the coefficients b. 
Interestingly, the ratio &mera/^mps for the four models produces very similar results, 
namely 1.80, 1.74, 1.72, and 1.56 for the Ising, Potts, XX and Heisenberg zig-zag models 
respectively. Given that the average ratio is &mera/&mps ~ 1-72, we conclude that in 
the large x limit a similar error in the energy for MERA and MPS is obtained if 

A^mera (x) « A^mps (X 1 - 72 ) , (43) 

that is, if Xmps = (Xmera) 1 ' 72 - Taking into account that the number of variational 
parameters in the MERA and MPS scale as (xmera) 4 and (xmps) 2 , this comparison 
shows that in the large x limit the MPS requires less variational parameters than the 
scale-invariant MERA in order to obtain a similar accuracy in the ground state energy. 
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TABLE II. Best fit coefficients to the functional form of Eq. [42] for the scaling of relative 
energy error in ground state MERA and MPS calculations. The central charge c of the 
critical models is given for reference. 





(i) MERA 
a b 


(ii) MPS 
a b 


Ising Model (c = 1/2) 
Potts Model (c = 4/5) 
XX Model (c = 1) 
Heisenberg Zig-Zag (c = 1) 


1.13 6.81 
17.80 5.22 
5.25 4.30 
1.89 3.80 


0.013 3.78 
0.432 3.00 
0.103 2.50 
0.059 2.44 



It is tempting to extend this comparison to computational costs. A first step in 
this direction is to note that each iteration in the optimization of MERA and MPS 
scales (naively) as (xmera) 6 and (xmps) 3 , from which it would be tempting to conclude 
that MPS algorithms require a lower computational budget than MERA algorithms to 
obtain the same accuracy in the ground state energy. However, there are important 
multiplicative prefactors fcMERA an d ^mps modifying the naive scaling of costs. In 
both cases, one is required to find the dominant eigenvector of a transfer matrix. But 
while in the case of the MERA this transfer matrix (or scaling superoperator) has a 
well-defined gap, implying that fcMERA is essentially independent of Xmera, in the case 
of the MPS the gap in the transfer matrix closes to zero with increasing Xmps i an d the 
prefactor fcMPS actually also grows with growing Xmps ■ Therefore a proper comparison 
of computational costs requires first a careful characterization of the dependence of 
fcMPS in Xmps, which is beyond the scope of the present manuscript. 

2. Two-Point Correlators 

Let us now compare the accuracy of two-point correlators produced by the scale- 
invariant MERA and MPS. For both approaches, we consider the ground state of the 
quantum XX model and compute the correlator, 

C (d) = (^(r)a(r + d)) (44) 

where a is a fermionic operator defined in terms of spin operators as 

a(r)=fllW) *V (45) 

\m<r / 

The correlation function of Eq. [44] has the exact expression 

C c ^ t (d)= - Sini f 2 \ (46) 

which is obtained from mapping the quantum XX model to a free-fermion model 4Bl ■ 
The correlation function C(d) decays polynomially, as expected from the ground state 
of a quantum critical model. 
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FIG. 12. (i) The two-point correlator of fermion operators, as denned Eq. |44| in the ground 
state of the quantum XX model, comparing results from the scale- invariant MERA and infinite 
MPS with the exact correlators. Correlators from MPS approximate polynomial decay only 
until some finite length scale £ « x ' 38 j while correlators from MERA remain polynomial at 
all length scales, (ii) Relative error in correlators as defined Eq. 48 



Fig. |T2|(i) shows the correlations obtained with a scale-invariant MERA for Xmera = 
{16, 46, 96} and with an MPS for Xmps = {24, 128, 512}. These particular values of the 
bond dimension \ have been chosen so that the tensor networks produce comparable 
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errors in the ground state energy, namely 

A-Emera (Xmera = 16) « ASmps (Xmps = 24) « 3.5 x 10~ 5 

A^mera (Xmera = 46) « AS M ps (Xmps = 128) « 4.5 x 10~ 7 

A^mera (Xmera = 96) ps ASmps (Xmps = 512) a 1.7 x 1(T 8 . (47) 

The figure illustrates quite clearly that while the MPS can only approximate polynomial 
correlations up to a finite length scale C = X K [33 , with k k, 1.38 for the quantum XX 
model, correlations in the MERA decay polynomially at all length scales [HI US] • Fig- 
[l2[ii) displays the relative error in correlators, 

A /T Inexact ^numeric I / ao\ 

AC - |7y 1 • I 48 J 

I ^ exact | 

Here it is seen that a MERA and an MPS that produce the same accuracy in the 
ground state energy produce also similarly accurate correlators at short distances, but 
the relative error in the correlator grows much slower with distance in the case of the 
MERA. For instance, the \ = 96 MERA reproduces correlators up to d = 10 6 sites 
with relative error AC < 2 x 10~ 4 , whereas the \ — 512 MPS, although possessing 
a similar ground energy, reproduces similarly accurate correlators only up to d f=a 100 
sites. 



3. Summary of Comparison 

To summarize, we have seen that the MPS is more efficient than the scale-invariant 
MERA when it comes to ground state energies of critical Hamiltonians, in that it re- 
quires less variational parameters to achieve the same accuracy. However, the MERA 
produces better correlators at large distances, and it is therefore better suited to char- 
acterize asymptotic behaviors, such as the polynomial decay of correlations, from which 
one could in principle extract the critical exponents of the theory. However, as exem- 



plified in Sect. IV D with the computation of scaling dimensions for local and non-local 
operators, critical exponents and other conformal data can actually be extracted more 
directly by analyzing the scaling superoperator. We illustrate this next. 



B. Evaluation of Conformal Data: The Ising Model 

As an example of extraction of conformal data from the scale-invariant MERA, here 
we identify the whole operator content (local and non-local primary fields) of the CFT 
corresponding to the quantum critical Ising model, reproducing the analysis of Refl2"2"l 
Similar results have also been obtained for the three-state Potts and quantum XX 
models in Refs. PHH2]. 

The Hamiltonian -ffi s i ng of Eq. [36|has a global internal Z2 corresponding to flipping 
all the spins. That is, Q = Z2 and~g G {+1, —1}, with V+i = I and V-\ = Z, and 



00 

r_l -fflsing r_j = iTlsingj T-l = Z. 

m— — 00 



(49) 
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The tensors u and w that comprise the scale-invariant MERA are chosen to be parity 
preserving; each index i of tensors u and w decomposes as i — (p, a p ), where p labels 
the parity (p = +1 for even parity and p = — 1 for odd parity) and a p labels the distinct 
values of i with parity p. For tensors u, w to be parity preservin g it is ensured that, 



if p{ii)p{i2)p{ji)p{ji) ~ ~ 1) m accordance with Eq. 16 An operator O 
acting on the spin chain has parity p if 

(r_oo(r_i)t=pO. (so) 

The local scaling superoperator <S g =i and the non-local scaling superoperator 5„ = _i, 



see Fig. 10'iv-v), are obtained from the optimized scale-invariant MERA with bond 
dimension \ — 32. The scaling superoperator <S g =i is diagonalized to find the local 
scaling operators 4>+i, a together with their scaling dimensions A + i, Q , while the scaling 
superoperator S g —-\ is diagonalized to find the non-local scaling operators of the form 



-X.a 



= ■ ■ ■ Z <g) Z ® Z (g) <f>-i, a , (51) 



together with their scaling dimensions A_i jQ . 

Table |Hl| compares the exact scaling dimensions of the primary fields of the Ising CFT 
and their numerical estimates obtained from the MERA, which reproduce the former 



with 4 to 6 digits of accuracy in all cases. In Fig. 13 we plot the scaling dimensions of 
magnitude A < 2.5 obtained from the scale-invariant MERA, which correspond to the 
primary fields and their descendants, organized both according to the locality g and 
the parity p of the corresponding scaling operators. Local scaling operators (g = +1) 
with even parity {p = +1) form the two conformal towers of the primary fields identity 
I and energy e of the Ising CFT, whereas those with odd parity (p = —1) form the 
conformal tower of the primary field spin a. Non-local scaling operators (g — —1) with 
even parity (p = +1) form the conformal tower of the disorder operator /i, and those 
with odd parity (p = —1) are organized according to two towers corresponding to the 
fcrmion operators ij) and ip. The numerical results from the scale- invariant MERA are 
seen to accurately reproduce the smallest scaling dimensions, those with A < 2.5, from 
all six conformal towers of the Ising CFT G2] • 



TABLE III. Scaling dimensions of the primary fields of the Ising CFT. 



A =xa=t 


A MERA 
= 32 


Error 


A CT =0.125 

A £ =l 

A M =0.125 

Av,=0.5 

A^=0.5 


0.1249998 
1.0001139 
0.1250002 
0.4999959 
0.4999963 


2 x 10~ 4 % 
0.011 % 
2 x 10" 4 % 
8 x 10" 4 % 
7 x 10- 4 % 



We have also computed the OPE coefficients C Q , / 3 7 for all primary fields, obtained 
by analyzing three-point correlators as described in Refs. [TB]. Table IV shows the 



numerical estimate of all non-vanishing OPE coefficients. Once again the results match 
their exact values to within several digits of accuracy. Thus, not only have we been 
able to identify the entire field content {I, e, er, ip, ip, fi} of the Ising CFT from a simple 
and rather inexpensive analysis of a quantum spin chain, but we can now also identify 
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FIG. 13. A few scaling dimensions of the critical Ising model obtained from a \ = 32 scale- 
invariant MERA. The scaling dimensions are organized by the both locality g = ±1 (local/non- 
local) and parity p = ±1 (even/odd) of the corresponding scaling operators and together form 
the six conformal towers of the Ising CFT. 



TABLE IV. OPE coefficients for the local and non-local primary fields of the Ising CFT. 





xact 


/— rMERA 
°X = 32 


error 




= 1/2 


0.50008 


0.016% 




= -1/2 


-0.49997 


0.006% 


c - 


B -i*/4 


1.00068e~ i7r/4 


0.068% 
0.068% 




1.00068e" r ' 4 


_ V2 






= i 


l.OOOli 


0.010% 


o - 


= —i 


-l.OOOli 


0.010% 



all possible subsets of primary fields that close a subalgebra by inspecting Table [TV] 
Indeed, it follows that we have the fusion rules 



€ x e 



I. 



a x a 



It- 



er x e 



/jx/i = i+£, (j, x e = (j., 

■0 X lp — I, tp X %i> = I, 

ip x ip = e, ip x e = ip, tp x e = ip, 



(52) 
(53) 
(54) 
(55) 



(as well as other, such asaXfi = ip + 'ip, etc) from where we see that {I, e} and 
{I, e, cr} close subalgebras of local primary fields, whereas {I, e, /i} and {I, e, tp, ip} close 
subalgebras that contain both local and non-local primary fields, where locality is 
relative to the spin variables. 
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VI. CONCLUSIONS 



In this manuscript we have presented an introduction to the scale-invariant MERA 
and its application to the study of quantum critical systems. The main strength of 
MERA, when applied to quantum critical systems, is that it can explicitly incorporate 
scale invariance. This facilitates enormously the computation of scaling dimensions 
(equivalently, critical exponents) and of other properties that characterize a quantum 
phase transition. 

Direct comparison with an MPS shows that, while the later is more efficient at 
computing local observables such as the ground state energy, the MERA produces 
significantly more accurate correlators at long distances. In addition, from the MERA 
it is straightforward to identify the scaling operators of the theory, as well as their 
scaling dimensions and operator product expansion, producing accurate conformal data 
that can be used to unambiguously identify the underlying CFT. 

Here we have considered homogeneous systems. However, the scale-invariant MERA 
has been successfully generalized to critical systems where translation invariance is 
explicitly broken, as it is the case of a critical system with a boundary, with an impurity, 
or the interface between two critical systems [^151] . In all these scenarios translation 
invariance is no longer present, but exploitation of scale invariance still produces a 
MERA algorithm with a cost 0(1) (that is, independent of the lattice size N), and 
therefore infinite systems can be addressed. Again, scaling operators associated to 
boundaries, defects and interfaces can be easily extracted from the simulations. 

Finally, much of the MERA formalism for critical systems in D = 1 dimensions is 
also directly applicable to D — 2 dimensions, including the characterization of scaling 
dimensions. However, due to significantly larger computational costs, so far only small 
values of x have been used in actual computations. Thus, further progress needs to 
be made in reducing computational costs before the scale-invariant MERA becomes a 
viable approach to quantum criticality also in D = 2 dimensions. 



Appendix A: Optimization algorithm for the scale-invariant MERA 

In this section we describe an algorithm to optimize the scale-invariant MERA to 
approximate the ground state of a critical system. It is based on modifying the opti- 
mization algorithm for a regular MERA of Ref. [TS] and has been previously sketched 
in Refs. [HI [16], although the present implementation differs in some details. 

A scale-invariant MERA is composed of a small number M of transitional layers 
{f/M , C/ [1] , . . . , U^ M ~^ } followed by an infinite sequence of (identical) scale-invariant 
layers, here denoted with an asterisk, U* . We shall henceforth use similar asterisk 
notation for all operators and tensors associated to the scale-invariant layers. Each of 
the transitional layers are characterized by a single isometry tuM and disentangler 
vy'] similarly, the scale-invariant layers U* are characterized by a single isometry w* 
and disentangler u* . The goal of this section is to describe how these tensors can be 
optimized in such a way that, given a critical Hamiltonian, the scale-invariant best 
approximates its ground state. 

We begin, in Sect. |A 1[ by describing the key building blocks of the optimization 
algorithm, which are (i) coarse-graining of the Hamiltonian; (ii) fine-graining of the 



density matrix; and (iii) optimization of one tensor of the MERA. Then, in Sect. A 2 
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we describe how these algorithmic building blocks can be put together to form an 



iterative optimization scheme. Finally, Sect. A3 describes computational tricks to 
improve convergence and accuracy. 





(iii)/^ Step.1: Initialize tensors 

Density Matrices: $p m ,fp^,p' 

Hamiltonian Terms: j y ff& a h™,}? } 

Disentanglers: Ui m ,U W ,u"\ 

Isometries: j^ 01 , W 1 ' 1 , W*} 



Step. 2: Update density 
matrices 



\p'",-^.p-]^\P",p>",p-] 



Step. 3: Update transition 
layers 



Step. 4: Update scale 
invariant layers 




FIG. 14. (i) Each layer C/' T ' of the MERA may be thought of as a coarse-graining transform 
between an initial lattice £' r ' and a coarser lattice £' T+1 '. (ii) A schematic representation 
of the scale-invariant MERA, here with two transitional layers C/' ' and followed by an 
infinite number of identical scale-invariant layers U* . (iii) The steps for optimizing a scale- 
invariant MERA as described Sect. IA 21 
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1. Building blocks 



a. Coarse-graining of the Hamiltonian 



A key part of the MERA optimization algorithm is the coarse-graining of the Hamil- 
tonian. Given an initial Hamiltonian = ~^2 r h^(r,r + 1) that decomposes as a 
sum of identical local terms /zi°](r, r + 1) = h° , we wish to construct, for any level r, 
the coarse-grained Hamiltonian = ^ r h^(r, r + 1), defined as 

H [r] = (ulr-l]V H [r-1] (tf[r-lrt . (A1) 

The coarse-graining of an operator that decomposes as a sum of local operators, such as 
a local Hamiltonian, can be achieved with the ascending superoperator formalism intro- 
duced in Sect. |II A| for the coarse-graining of a single local operator. The coarse-grained 
Hamiltonian coupling fti r+1 J is obtained by enacting the (left, center, right) ascending 
superoperators Al,Ac and Ar on the Hamiltonian term W- T \ 

= A l [ ] (h^) +AP (ftM) +4 Tl (hW) 

= 3_4 M (V rl ) (A2) 

where A is the average of the three ascending superoperators. The diagrammatic 



representation of the tensor network described by Eq. A2 is shown in Fig. 15 'i) 



b. Fine-graining of the density matrix 



Another key part of the MERA optimization algorithm is fine-graining the density 
matrix. Given the average two-site density matrix defined on lattice C^ T \ we 
wish to construct the average two-site density matrix pt T_1 l on lattice £[ T_1 1. This 
is accomplished through the average descending superoperator T> introduced in Sect. 
|IV C[ which acts as 



p[r-l] = ^[r-l] 



(A3) 



From iteration of Eq. 



A3 



for all 
Eq 



one can construct the average two-site density matrices p' T > 
The diagrammatic representation of the tensor network described by 
A3 is shown in Fig. |15|(ii). 



T < T. 



c. Optimization of one tensor of the MERA 



In order to approximate the ground state of a Hamiltonian H 7 the disentanglers u 
and isometries w that define a scale- invariant MERA |W) should be chosen to minimize 
the energy E = (^\H\^>) . We shall proceed by updating one tensor of the MERA at 
a time, while holding the rest of the tensors fixed. Given a tensor to be updated, we 
note that the energy E depends quadratically on that tensor and its conjugate. Here 
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(i) Lifting the Hamiltonian 




I* ♦ 



i r 



1 O l 

w 




-vr 



(ii) Lowering the Density Matrix: 




(iii) Environment of an isornetry: 

T 




T [«-l] 








< 

















(iv) Environment of a disentangler: 

Y 






FIG. 15. This figure displays the full set of tensor network contractions required to optimize 
a ternary MERA. (i) The tensor network contractions required to coarse-grain a local Hamil- 
tonian, see Sect. Ala (ii) The tensor network contractions required to fine-grain the average 
two-site density matrix, see Sect. Alb (iii) The linearized environme nt T w of an isornetry 
w and (iv) the linearized environment T u of a disentangler u, see Sect. Ale 
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we employ the linearized optimization scheme of Ref. |15j . which we sketch next (a 
justification for linearizing the cost function and further details can be found in Ref. 

USD- . 

The update of an isometry w relies on computing its environment T w , which repre- 
sents a factorization of E = (^\H\ 1 ii) to isolate its dependence on w, 

E = tr(wT w ) + k u (A4) 

with ki an irrelevant constant. [T w is the derivative of (^IH^) with respect to w 
while keeping fixed] . The updated isometry w' that minimizes the linearized energy 
— that is, such that tr(w^ T m [ T ]) is minimal — is given by w' = — V^V^, where V\ 
and Vi are obtained from the singular value decomposition (SVD) of the environment, 
T w — V1SV2 ■ We then proceed by replacing that particular isometry w with w' in 
the MERA, which now represents a new state j*'). We emphasize that the energy 
E' = (^'\H\^>') is computed by replacing both w and uft with w' and w'\ and it 
is therefore not given by tr(w'Y w ) + k\. In other words, w' does not minimize the 
expectation value of the Hamiltonian (in fact, the energy could even rise!). However, the 
approach works well in practice. The update of a disentangler u follows an analogous 
procedure, involving the SVD of the corresponding environment T u . 

The environments of an isometry w^ T ' and disentangler at layer r depend only 
on a small number of other tensors, 

T roM = T wM (uH «;M pt T+1 U W ) . (A5) 
T„ w =T uW ( M M W [T] ,P [T+1 U M ), (A6) 



as depicted in Figs. 15 iii)-(iv) 



2. Optimization algorithm 

For concreteness, we describe the optimization of a scale-invariant MERA with M = 
2 transitional layers as depicted Fig. [Ti^ii). The ansatz is then completely characterized 
by three isometries , w* } and three disentanglers {u^, u* } , where the 

isometry w* and the disentangler u* define the scale-invariant layers. For purposes 
of the optimization algorithm, one should also store in computer memory the two-site 
Hamiltonian couplings {h}°\ ftW, hffl, h^\. Here is a scale-averaged Hamiltonian, 
to be introduced later in Eq. |A11| We also store the average two-site density matrices 
{pl 1 !,^*}, where we recall that p* is the two-site density matrix any of the scale- 
invariant layers. We now proceed to describe the optimization algorithm for scale- 



invariant MERA, as outlined Fig. 14 'iii) . In an iteration of the algorithm, comprising 
Sect. |A2b| Sect. |A 2 c| and Sect. |A2d[ each of the stored tensors is updated to a new 
version of itself denoted with a prime, e.g. 1— > iJ !' . 



a. Initialization 



The preliminary step of the algorithm is to initialize the tensors which define the 
MERA. It is most often sufficient to initialize the tensors to be random, though in 
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certain situations, for instance if one has prior knowledge about the ground state of the 
Hamiltonian, it is useful to use that information. Random isometries {id' ', ly' 1 ', w*\ 
can be obtained through singular value decomposition of an appropriately sized random 
rectangular matrix. The disentanglers {it' ', it' 1 ' 1 u* } can also be initialized randomly, 
or as the identity operator. Notice that the local Hamiltonian term ft' ' is given as 
the input of the algorithm (it describes the critical Hamiltonian whose ground state 
is to be found), whereas the Hamiltonian terms {ft.' 1 ', ft.' 2 ', ft^} and density matrices 
{p' 1 ',^*} will be generated during the first iteration of the algorithm. 



b. Updating the density matrices 



The first step of each iteration is to update the two-site reduced density matrices. 
As explained Sect. IV C the scaling density matrix p*' is defined as the dominant 
eigenoperator of the average descending superoperator f>* associated to the scale- 
invariant layers, 



P 



V* (p* 



(A7) 



One can solve Eq. A7 for p*' with a sparse diagonalization technique such as the 
Lanczos method, where the scaling density matrix p* from the previous iteration (if 
any) can be reused as the starting point for the Lanczos method to accelerate conver- 
gence. The updated density matrix p' 1 '' is then obtained by fine-graining p*' with the 
descending superoperator 2?' 1 ', 



pM' = £[l] (p*'). 



(A8) 



c. Updating the transitional layers 



The next step of the iteration is to update the transitional layers of the MERA. The 
tensors u>' ' and it' ', which comprise the first transitional layer t/' ', are updated by 
first computing their linearized environments, 

r ul0] =r um (uM,wM,pW,hM), (A 9 ) 

from which the updated tensors ui' '' and it' '' are obtained through SVD, see Sect. 



Ale Together the updated isometry it;' ' ' and disentangler it' '' define the ascending 
superoperator ^4'°', which is used to coarse-grain the Hamiltonian ft' ', 



ft' 1 '' = 3A [0] (ft' 01 ) , (A10) 



and obtain the updated Hamiltonian coupling ft' 1 '', as described Sect. Ala The 
second transitional layer t/' 1 ' is then updated, in the same manner as described here 
for the first transitional layer, to obtain updated tensors {it;' 1 '', it' 1 '', ft' 2 ''}. 
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d. Updating the scale-invariant layers 



The final step of the iteration is to update the isometry w* and disentangler u* asso- 
ciated to the scale-invariant layers of the MERA. Updating the scale-invariant layers is 
more complicated than updating the transitional layers since the scaling tensors w* and 
u* do not appear in just one layer but in all layers for r > M, with M the number 
of transitional layers). To update these tensors, we construct linearized environments 
that average the contributions coming from all length scales r > M. Computation of 
these scale-averaged environments is simplified [15] by the introduction of the scale- 
averaged Hamiltonian term , which for the case of M transitional layers is defined 
as 

h<> = £ hM = hM + \ hW+l] + l hlM+2] + ■ ■ ■ > (All) 

T = M 



where the factor of (1/3) r M arises from the fact that layer [/M contains three times 
as many tensors as layer U^ T+1 ^ |15j . 

Note that for a critical Hamiltonian (with negligible or sufficiently suppressed RG- 
irrelevant terms) and a scale invariant MERA that has already been optimized, scale- 
invariance implicitly assumes that /)i T+1 J oc for all r > M, and therefore the 
scale-averaged Hamiltonian is simply proportional to fti M l . However this property 
will only hold for a MERA which has been properly optimized; during the optimiza- 
tion procedure it is necessary to explicitly average over all scale-invariant layers when 
constructing environments. Computing the scale-averaged Hamiltonian directly 



through Eq. All is not feasible due to the infinite summation. O ne str ategy to obtain 



an approximation to hfc is to use a partial summation of Eq. All which stops at 
some t = T. In practice, however, a useful estimate of the updated scale-averaged 
Hamiltonian hr' is already obtained from the operator from the previous iteration 
through 

hf>' w ftM' + A* (h?) . (A12) 

This estimate of hr' is accurate provided that the effective Hamiltonians are only 
changing by small amounts between iterations, i.e. h^' ~ for all t, which becomes 
a better approximation as the optimization nears convergence. Once the new h^' has 
been computed then the environments of the scaling isometry w* and disentangler u* 
are computed, 

TV = r w * (u*, w *,p* , ,h^ / ), 

T u ,=T u ,(u*,w*,p* , ,hO , ), (A13) 



from which, as usual, one obtains the updated tensors w*' and u*' by SVD. This 
concludes a single iteration of the optimization algorithm for scale-invariant MERA; 
the next iteration can begin again at Sect. |A 2 b| 
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FIG. 16. (i) The preliminary blocking groups two d-dimensional sites of the original lattice 
into a single d 2 -dimensional site of the lattice £' '. (ii) The tensor it* 1 ; is equiva- 
lent to the identity, see Eq. A14 (iii) Under the preliminary blocking the initial Hamil- 
tonian = 5~J r ft' ' (r, r + 1) on lattice £' ' is mapped to an equivalent Hamiltonian 

H l0] = E r h [0] (r, r + 1) on lattice £ [0] . 



3. Computational tricks 



a. Preliminary blocking 



Given a local Hamiltonian defined on a lattice , where each lattice site is described 
by a Hilbert space V of finite dimension d, in certain cases it may be convenient to 
perform a preliminary blocking of e.g. two sites of into a single site of dimension 
d = d 2 in the coarser lattice CM. This preliminary blocking, which maps the initial 
Hamiltonian =Y] on lattice to a new Hamiltonian = ^ defined 
on a coarser lattice C<-°>, simply amounts to re-expressing the Hamiltonian in a different 
form, i.e. and H<- ' are different representations of the same Hamiltonian. An 

example how the Hamiltonian can be redefined through a preliminary blocking, for 
the specific case of blocking two sites into one, is depicted Fig. 
equivalent to the identity, 



16 with u a tensor 




(A14) 



where indices %2 and «3 take d values each whereas index i\ takes d 2 values. This 
preliminary blocking, which does not increase (to leading order) the contraction cost 
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of the MERA provided d 2 is less than the bond dimension \ used in higher layers of the 
MERA, has two potential advantages. Firstly, it can reduce a next-nearest neighbor 
Hamiltonian into a nearest neighbor Hamiltonian (of larger local dimension), 
which can then be easily treated with the ternary MERA. The second advantage is 
that it transforms a state that would otherwise be translation-invariant by shifts of 
two sites into a state that is translation-invariant by shifts of just one site, which can 
then be represented directly with a translation-invariant ternary MERA. 



b. Shifting the Hamiltonian spectrum 



Let us assume we are interested in obtaining the ground state MERA of a given 
local Hamiltonian = ~J2 r h\°\r, r + 1). Prior to the optimization it can be useful 
to shift the spectrum of the Hamiltonian by adding or subtracting contributions of the 
identity, 

i— > h® = fcP>] — al, (A15) 
where I is the identity operator on two sites, such that the shifted Hamilto nian H$ = 



^ ha is negative defined. This can be achieved by choosing a in Eq. 



A15 



as the 



largest eigenvalue of /J '. Having a negative defined Hamiltonian ensures that the 
optimization targets the low-energy, i.e. ground state, subspace (the optimization 
algorithm, based upon extremizing the energy of the state, could otherwise target the 
high-energy subspace). 

It is also useful to similarly shift the spectrum of the coarse-grained Hamiltonian 
that arise during the optimization, by replacing i— >• where a as the largest 
eigenvalue of h\' T \ since this is seen to speed-up convergence during the optimization 
of the MERA. 



c. Converging the number of transitional layers 



In order to obtain an accurate approximation to the ground state of critical Hamil- 
tonian H with a scale invariant MERA it is important to use enough transitional layers 
to sufficiently suppress the effect of any RG irrelevant terms present in H, which can 
break scale- invariance at short length scales. For a given critical H the appropriate 
number transitional layers is not known a priori, and must be found by trial and error. 
We proceed in the following way. Suppose that we have already optimized a MERA 
with M = 2 transitional layers, that is, characterized by {U^°\ U*}. We can 
then add an additional transitional layer, M = 3, so that now the MERA is given by 
{/J! ', C/W, U^ 2 \ U*\, where we initially set = U* , and re-optimize the tensors as 
to best approximate the ground state of Hamiltonian H. We then keep adding more 
transitional layers, and re-optimizing the ansatz, until addition of an extra transitional 
layer does not produce significantly different results under re-optimization in terms e.g. 
of the ground energy, scaling dimensions, etc. 
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Appendix B: Comparison of MERA schemes 



In Sect. II C three different implementations of the MERA in D = 1 dimensions were 
described: the binary MERA, the ternary MERA and the modified binary MERA. The 
three schemes differ in how the computational cost (incurred e.g. in optimizing the 
ansatz and in computing the expectation value of local observables) scales with the 
bond dimension %, namely as 0(x 9 ), 0(x 8 ) and 0(x 7 ), respectively. In addition, the 
three schemes also differ in the strength of disentangling, which means that for the same 
value of x the numerical accuracy of the results will vary across the three schemes. In 
this section we investigate which of the schemes provides the most accurate ground 
state energy for a fixed computation cost. 

For this purpose, we have optimized each of the three MERA schemes so as to 
approximate the ground states of the critical Hamiltonians defined Eqs. [36fl39| in Sect. 



|V) Fig. 17 shows the relative energy error AE, defined as 

AE = (-Eexact -^numeric )/-^exact 3 (-^-0 

as a function of bond dimension x- For all three schemes the relative energy error AE 
is seen to scale as a power of the bond dimension x, 

AE « a X -\ (B2) 

for some coefficients 'a' and '6' that depend both on the MERA scheme in use and 
the spin model under consideration. The sets of these coefficients, obtained from 
linear fits in Fig. [T7J are displayed in Table [V] For each of the four spin models 
separately, the data shows that while the coefficient 'a' depends considerably on the 
MERA scheme (with variations by one or two orders of magnitude), the coefficient 
'6' is remarkably constant across the schemes (with variations of about 15%). These 
results indicate that, while using a scheme with more powerful disentangling (such 
as the binary MERA scheme) can reduce the energy error by a considerable, fixed 
multiplicative factor (as expressed by a significantly smaller constant 'a'), the power 
'6' of x f° r a given critical model is not significantly affected by the disentangling power 
of the scheme, even though the latter does affect significantly how the computational 
cost scales with x- This suggests that the modified binary MERA scheme, which has 
the smallest disentangling power, will be the most cost effective at large Xi since its 
cost scales as a smaller power of x an( i the accuracy scales roughly as the same power 
of x as the other schemes. 

Fig. fl8] shows AE directly as a function of the computational cost for the Ising model. 
Similar results are obtained for the other models. The costs have been approximated 
to be C = x 9 for binary, C = x 8 for ternary and C = x 7 for modified binary MERA, 
where C is measured in floating-point operations. The exact expression for the cost 
would include a multiplicative constant of order one that does not change significantly 
across the schemes and has been ignored for simplicity. The plot confirms that, for 
large Xi the modified binary MERA scheme, with cost C — x 7 , is the scheme that 
provides the most accurate results for a fixed computational cost, whereas for smaller 
values of x the binary scheme is the most cost effective. The crossover occurs at a 
cost of order C ~ 10 7 floating-point operations per iteration, which is well within the 
capabilities of a small workstation. Therefore, in most practical numerical applications 
the modified binary MERA is a better scheme than the others. 
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FIG. 17. Relative energy error AE in the ground state of critical spin chains, as a function 
of the MERA bond dimension x, comparing three different MERA schemes. 



TABLE V. Best fit coefficients to the functional form of Eq. |B2| for the scaling of relative en- 
ergy error in critical ground states, comparing the three different MERA schemes as described 
in Sect. |II C| The central charge c of the critical models is given for reference. 





(i) binary 
a b 


(ii) ternary 
a b 


(iii) mod. bin. 
a b 


Ising Model (c = 1/2) 
Potts Model (c = 4/5) 
XX Model (c = 1) 
Heisenberg Zig-Zag (c = 1) 


0.012 6.20 
0.527 4.63 
0.615 3.97 
0.824 3.90 


0.105 5.98 
1.786 4.58 
1.156 4.05 
3.264 3.95 


1.13 6.81 
17.80 5.22 
5.25 4.30 
1.89 3.80 
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FIG. 18. Relative energy error scaling in ground state MERA calculations, for the three 
different MERA schemes as described in Sect. |II C| in the quantum critical Ising model 
plotted as a function of the leading order computational cost C for the optimization. 



Appendix C: Reducing the Cost of MERA Contractions 

In Appendix [Bj a numerical study determined the modified binary MERA scheme, 
which can be optimized with a leading computational cost which scales as 0(x 7 ), to 
be the optimal ID MERA implementation of those considered. In this Appendix we 
describe how the cost of implementing the modified binary MERA scheme can be 
reduced from 0(x ) to 0(x 6 ) through the use of approximations in the tensor network 
contractions required for its optimization. The steps required to optimize the modified 
binary MERA are analogous to those described for the ternary MERA in Sect. [X] 



Fig. 19 shows four closed tensor networks (tensor networks without any open index). 
The tensor networks required to optimize the modified binary MERA (which include 
the ascending and descending superoperators, and the environments for single tensors) 
can be generated from these closed tensor networ ks t hrough removal of a single tensor. 



The cost of contracting the two networks in Figs. 19 'i-ii) scales as 0(x 7 ), whereas cost 
of contracting the two networks in Figs. 19 iii-iv) scales just as 0(x 6 ) contraction cost. 



If we could reduce the cost of the two first networks also down to 0(x 6 ), then that 
would be the overall leading cost of the modified binary MERA. Let us then see how 
this can be accomplished. 



1. Insertion of Projectors 



Let us consider a rank x projector P, decomposed as the product of an isometric 
tensor v and its conjugate v* , 

P = W,\ ^.(VJ^^ (C1) 
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FIG. 19. (i-iv) The four types of tensor network required for the optimization of the modified 
binary MERA. In terms of bond dimension x> the networks of (i) and (ii) are contraction 
cost 0(x 7 ), while networks of (hi) and (iv) are of contraction cost 0(x 6 )- (v-vi) The tensor 
networks of (i) and (ii) have been modified with the inclusion of a rank x projector P = 
vv^ . The modified networks are of contraction cost 0(xx 5 )- ( vu ) The tensor «' T ' should be 
optimized such that projects onto the subspace of the density matrix p' T ' with 

greatest weight, as described Eq. |C2| 



where V x is a ^-dimensional vector space and is a x-dimensional vector space 
for some x ^ X 2 ■ I n place of using the original tensor networks of Fig. 19 i-ii) 



in the optimization algorithm, with contraction cost 0(x 7 ), we shall use the tensor 
networks of Fig. 19 V-vi), which have been modified through inclusion of the projector 
P and whose cost is now of order 0(x 5 X)- I n the reminder of this Appendix we first 
explain how to optimally choose the isometric v so that the modified tensors networks 
yield environments that best approximate the environments obtained from the original 
networks; and then we argue, based on numerical evidence, that x mav indeed be 
chosen as O(x) without significant loss of accuracy, thus reducing the overall cost of 
optimization down to 0(x 6 )- 
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FIG. 20. The spectrums of one-site and two-site density matrices, obtained from the scale- 
invariant fixed point of a \ = 16 MERA optimized for the ground state of the quantum XX 
model. The two site density matrix can be truncated down to its \ = 50 most significant 
eigenvalues, out of the total \ 2 ~ 256 eigenvalues, without significant loss of accuracy (as 
gauged by the smallest eigenvalue of the one-site density matrix). 



2. Optimizing the Projector 



Given the two-body density matrix pM at level of the MERA, the proper choice 
for the isometric is that which maximizes the trace of the density matrix, 



lax ( tr ( 

,m V V 



]t p M w M 



(C2) 



as we now justify. Under this choice of the projector = yMuMt projects 
onto the subspace of pM spanned by its eigenvectors of largest eigenvalue. If we then 
assume that only \ eigenvalues out of the total x 2 eigenvalues of pM have significant 
weight (the rest being negligibly small) then the density matrix is invariant under the 
projection, pM = pMpM, to within, by assumption, negligible error. Thus, in this 
case, when computing the density matrix from the higher level density matrix 
p[r+i] ^ modified tensor networks of Fig. 
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i) will give the same result as the 
original tensor networks of Fig. 19 i-ii). Likewise the other tensor environments (such 



as environments of disentanglers and isometries) generated from the modified tensor 
networks will be the same as those generated from the original networks. 

In principle the iso met ric u' r ' that projects onto the most significant subspace of 
p^ , as described Eq. C2 could be obtained directly from the spectral decomposition 
of the density matrix p^ 1 " 1 for each t; however, computing the exact reduced density 
matrix pM from p[ r+1 l is an 0(x 7 ) operation and should be avoided. By expressi ng t he 



density matrix pI T l in terms of p^ T+1 \ and w^ T \ we can write the trace in Eq. C2 
the tensor network in the left hand side of Fig. 19 vii). From that tensor network, we 
can extract the environment for v with cost 0(x a)- Therefore, by using the linearized 



single tensor optimization described in Sect. |A 1 c we can iteratively optimize v at a 
cost 0(x 5 x)- 
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3. Rank of Projector 



The introduction of a properly optimized rank \ projector P was argued to be 
equivalent to truncation of the two-site density matrix to retain its x most significant 
eigenvalues; this leads to truncation error e = 1 — tr (pMpM) that should be kept 
sufficiently small. A good indication of exactly how small e should be kept is to ensure 
that it is the same size, or smaller than, the smallest eigenvalue of the one-site density 
matrix, which is indicative of the degree of accuracy of the MERA under consideration. 
Thus the relevant question becomes, for a MERA with bond dimension x, what is the 
necessary rank \ of the projector P required to maintain this sufficient degree of 
accuracy? 

Since the MERA represents states with at most a logarithmic scaling of the entropy 
with block size, see Sect. IV A the entropy of the two-site density matrix is much less 
than twice the entropy of the one-site density matrix. Therefore, it is to be expected 
that x < X 2 - [Notice that x — X 2 would be consistent with an entropy that scales 
linearly with the block size]. Indeed, the numerical evidence suggests that x may be 
chosen as 0(x)\ for the simulations with the modified binary MERA scheme of Sect. 
|Vj it was found that keeping x ~ 3x gave a sufficient level of accuracy for all four 
critical models under consideration, and over the large range of bond dimensions x 
analyzed. Given that this relation held over a range of bond dimensions between x = 4 
and x = 150, it seems likely that the relation between x an d X is indeed linear, or at 
least very close to linear for ground states of critical systems. An example is shown in 
Fig. [20j which plots the spectrum of one-site and two-site density matrices for x = 16 
quantum XX model, where it can seen that choosing x ~ 3x is sufficient to ensure that 
the truncation error e is of the same order as the smallest eigenvalue of the one-site 
density matrix. 

Thus, though not rigorously justified, the available evidence indicates that x may 
be chosen as 0(x), hence the overall cost of optimizing the modified binary MERA 
has been reduced to 0(x 6 ) cost. The reduction in cost comes at the price of both 
introducing new tensors that must be updated with each iteration and introducing 
a controlled amount of approximation into the tensor network contractions. Although 
we have only described how the cost of the modified binary MERA scheme can be 
reduced through introduction of a projector P into the tensor network diagrams, the 
same approach could be employed to potentially reduce the cost of any MERA scheme. 
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